class: inverse,left, middle background-image: url(data:image/png;base64,#background.png) background-size: cover <img src="data:image/png;base64,#LOGO_DIPLOMADO.png" width="500px"/> ##Módulo 2: EstadÃstica Espacial y geoestadÃstica ###GeoestadÃstica: Análisis Exploratorio de Datos Espaciales José A. Lastra<br> <a href="http://github.com/JoseLastra"> Github: JoseLastra</a><br> <a href="mailto:jose.lastra@pucv.cl"> jose.lastra@pucv.cl</a><br> .large[<b><a href="https://www.pucv.cl/uuaa/site/edic/base/port/labgrs.html">LabGRS</a> | Octubre, 2023 </b>] <br> --- class: center,middle background-image: url(data:image/png;base64,#labgrs_logo.png) background-size: 35% --- ## Contenidos .pull-left[ 1) GeoestadÃstica: Kriging 2) Variograma 3) Etapas del modelado geoestadÃstico: * Análisis Exploratorio de Datos Espaciales (AEDE) * Análisis estructural * Predicciones * Validación de modelos ] .pull-right[ <img src="data:image/png;base64,#https://raw.githubusercontent.com/allisonhorst/stats-illustrations/main/rstats-artwork/r_rollercoaster.png" width="650px"/> ] --- ## GeoestadÃstica: Kriging -- - La **geoestadÃstica** (Matheron 1962), presenta un conjunto de herramientas útiles para analizar y predecir valores de una variable en el dominio espacial o espacio-temporal. -- - Pertenece a la disciplina de la estadÃstica espacial, permitiéndonos describir y predecir la continuidad espacial de los datos de entrada. -- - Además de entregarnos la predicción de los valores, nos permiten obtener medidas de incertidumbre asociadas a la predicción. -- - Dentro de los métodos geoestadÃsticos **Kriging** es uno de los más empleados. <center><img src="data:image/png;base64,#https://gisgeography.com/wp-content/uploads/2017/01/kriging-prediction-standard-error-white-850x390.png" width="500px"/></center> <center><a href="https://gisgeography.com/">GISGeography, 2022</a></center> --- ## GeoestadÃstica: Kriging -- - La ventaja de estos sobre otros métodos consiste en que considera la **tendencia espacial** (variación de gran escala) y la **autocorrelación espacial** (variabilidad de pequeña escala). Esto expresado normalmente como: $$ `\begin{equation} Z(\mathbf{s}) = \mu + \varepsilon(\mathbf{s}), \end{equation}` $$ Donde: - `\(Z(\mathbf{s})\)` es el valor en una posición `\(s\)` - `\(\mu\)` representa la tendencia (promedio global) - `\(\varepsilon(\mathbf{s})\)` proceso de error de media cero que incorpora la dependencia espacial (residual). --- ## GeoestadÃstica: Kriging -- - Hay diferentes modelos de kriging con distintos supuestos acerca de `\(\mu\)` [(Fernández y Cotos, 2021)](https://rubenfcasal.github.io/estadistica_espacial/): 1. **Kriging simple (KS)**: Asume una media conocida y asume que el variograma existe y es conocido. 2. **Kriging ordinario (KO)**: Supone una media constante pero desconocida. Asume también, que por lo menos existe el variograma y este es conocido. 3. **Kriging Universal (KU)**: Supone una media desconocida y no constante en el espacio. Sin embargo, asume que esta es una combinación lineal (desconocida) de `\(p + 1\)` variables explicativas. También asume un variograma eistente y conocido. <center><img src="data:image/png;base64,#krigings.png" width="1300px"/></center> <center>Adaptado de: ESRI, 2021, Fernández y Cotos, 2021</center> --- ## Variograma -- - La clave de la modelación con kriging radica en la construcción del **variograma**, que es calculado para todos los pares de puntos serparados por una distancia (**h**). -- - El variograma es otra forma de ver la autocorrelación espacial de los pares de puntos según la distancia (Cañada, 2006, Fuenzalida et. al., 2015). .pull-left[ <center><img src="data:image/png;base64,#variogram.png" width="600px"/></center> <center>LABGRS, 2022</center> ] .pull-right[ - **Nugget**: Error aleatorio asociado al muestreo. Representa una discontinuidad puntual del semivariograma en el origen. - **Sill (Meseta)**: Zona del semivariograma donde la variación se nivela. - **Rango**: distancia en la que las observaciones son independientes ] --- ## Etapas del modelado geoestadÃstico -- - La construcción de nuestro variograma (**experimental**) nos permitirá luego ajustar un modelo para la creación de las predicciones. -- - El modelado del semivariograma es un paso crucial entre la descripción espacial y la predicción espacial. -- - En general, un proceso ideal para la modelación geoestadÃstica empleando Kriging organizarse de la siguiente manera: 1. **Análisis Exploratorio de Datos Espaciales (AEDE)**: consiste en explorar y preparar los datos previo a la modelación espacial y construcción del variograma. 2. **Análisis estructural**: consiste en el estudio de la continuidad espacial de la variable. En este paso, es cuando se realiza la construcción del **variograma experimental** y el ajuste del variograma teórico seleccionado. 3. **Predicciones**: Creación de las estimaciones de la variable en los puntos no muestrales. Esto considera la estructura de correlación espacial seleccionada (**variograma teórico**). 4. **Validación**: Evaluación de los errores de predicción obtenido. -- - En esta clase haremos una revisión del primero punto: **AEDE** --- ## Etapas de la modelación: AEDE -- - Para proceder a la modelación con kriging ordinario, tenemos que revisar si nuestros datos cumplen con algunos criterios: * Distribución normal de los datos (*) * Datos deben ser estacionarios o pseudo-estacionarios * Sin tendencia -- - Para hacer este ejercicio trabajaremos con la tabla de nombre **temp_promedio_cauquenes.csv** ```r #lectura temp_cauquenes <- read_csv('temp_promedio_cauquenes.csv') #paso a SF temp_cauquenes_sf <- st_as_sf(temp_cauquenes, #tabla coords = c('X','Y'), #coordenadas crs = 32719) #sistema de referencia ``` ---
--- ## AEDE: Normalidad de los datos -- - La mayorÃa de los métodos geoestadÃsticos (no todos) son óptimos cuando los datos siguen la **distribución normal** (forma de campana de Gauss) y KO no es la excepción. <center><img src="data:image/png;base64,#https://www.investopedia.com/thmb/WzHphax5NQFPMmdgDACQNs9s8LM=/1500x0/filters:no_upscale():max_bytes(150000):strip_icc():format(webp)/dotdash_final_Optimize_Your_Portfolio_Using_Normal_Distribution_Jan_2021-04-a92fef9458844ea0889ea7db57bc0adb.jpg" width="900px"/></center> <center>Sabrina Jiang © Investopedia 2021</center> --- ## AEDE: Normalidad de los datos -- - **Importante**: en este momento **NO** estamos revisando las caracterÃsticas espaciales. ```r ## histograma de los datos hist(temp_cauquenes_sf$t_mean1, xlab = 'LST °C medido') ## qqplot ggpubr::ggqqplot(temp_cauquenes_sf$t_mean1,ylab ='LST °C medido') ``` .pull-left[ <img src="data:image/png;base64,#DIPGEOPR_02_5_files/figure-html/unnamed-chunk-5-1.png" width="100%" /> ] .pull-right[ <img src="data:image/png;base64,#DIPGEOPR_02_5_files/figure-html/unnamed-chunk-6-1.png" width="100%" /> ] --- ## AEDE: Normalidad de los datos -- - Otras formas sencillas de evaluar normalidad son (Fuenzalida et al, 2015): * Coeficiente de Sesgo (skewness, CS, Webster y Oliver, 2001, Buzai y Montes, 2021): + 0 < |CS| < 0.5 Se acepta la función de probabilidad como normal + 0.5 < |CS| < 1 Es necesario transformar datos para normalizar; aplicación de raÃz cuadrada + |CS| > 1 Transformación de tipo logarÃtmica * Kurtosis: si los valores se mueven entre -3 y 3 se acepta la normalidad. .pull-left[ <center><img src="data:image/png;base64,#https://revistas.tec.ac.cr/index.php/eagronegocios/article/download/4456/4046/12439" width="400px"/></center> <center>Aguilar, 2019</center> ] .pull-right[ <center><img src="data:image/png;base64,#https://www.allaboutcircuits.com/uploads/articles/understanding-the-normal-distribution-parametric-tests-skewness-and-kurtosis-rk-aac-image2.jpg" width="400px"/></center> <center><a href="https://www.allaboutcircuits.com/technical-articles/understanding-the-normal-distribution-parametric-tests-skewness-and-kurtosis/">Keim, 2020</a></center> ] --- ## AEDE: Normalidad de los datos -- - La forma más directa es emplear algún test de normalidad, muchos de los cuales están disponibles en R base. -- - La hipótesis nula de la mayorÃa de estos test es que la **distribución de la muestra es normal**. Por lo que, si ***p < 0.01*** se puede decir que los datos **no siguen una distribución normal**. ```r # shapiro-Wilk test R base considerando p < 0.01 shapiro.test(temp_cauquenes_sf$t_mean1) ``` ``` ## ## Shapiro-Wilk normality test ## ## data: temp_cauquenes_sf$t_mean1 ## W = 0.97656, p-value = 0.01146 ``` --- ## AEDE: Normalidad de los datos -- - También disponemos de otras herramientas de librerÃas como [**moments**](https://CRAN.R-project.org/package=moments) o [**nortest**](https://CRAN.R-project.org/package=nortest) ```r # Jarque-Bera test considerando p < 0.01 moments::jarque.test(temp_cauquenes_sf$t_mean1) ``` ``` ## ## Jarque-Bera Normality Test ## ## data: temp_cauquenes_sf$t_mean1 ## JB = 3.9031, p-value = 0.1421 ## alternative hypothesis: greater ``` ```r # Kolmogorov-Smirnov test considerando p < 0.01 nortest::lillie.test(temp_cauquenes_sf$t_mean1) ``` ``` ## ## Lilliefors (Kolmogorov-Smirnov) normality test ## ## data: temp_cauquenes_sf$t_mean1 ## D = 0.077545, p-value = 0.02778 ``` --- ## Estacionariedad y autocorrelación espacial -- - Uno de los supuestos base del kriging ordinario es la **estacionariedad espacial débil** o **estacionariedad de segundo orden** -- - Esto definido básicamente implica que el promedio y la varianza son constantes (en una vecindad) y que la autocorrelación depende solo de la distancia entre las unidades muestrales. -- - Varias cosas pueden vulnerar el concepto de estacionariedad: * No estacionariedad de la media * No estacionariedad de la varianza * AnisotropÃa -- - **No existe test para evaluar la estacionariedad**, aunque se pueden crear visualizaciones como la de polÃgonos de Voronoi, variogramas direccionales, análisis de correlación a diferentes distancias u otros. -- - **Si podemos evaluar la autocorrelación espacial**, ya sea a nivel global o local (Moran I, variogramas, etc.) --- ## Estacionariedad y autocorrelación espacial -- - Primero evaluaremos si nuestros datos presentan autocorrelación espacial global empleando Moran I. -- - Calcularemos nuestros vecinos empleando bandas de distancia donde al menos haya un vecino para cada punto (umbral crÃtico) -- ```r temp_cauquenes_sf <- temp_cauquenes_sf %>% mutate( nb = st_dist_band(geometry), wt = st_weights(nb) ) ``` -- - Podemos emplear otras formas de obtener los vecinos como **st_knn()** o **st_contiguity()**. Cada una referida a puntos y polÃgonos. --- <img src="data:image/png;base64,#DIPGEOPR_02_5_files/figure-html/unnamed-chunk-10-1.png" width="90%" style="display: block; margin: auto;" /> --- ## Estacionariedad y autocorrelación espacial -- - Cálculo de global Moran I ```r nb <- temp_cauquenes_sf$nb #vecinos wt <- temp_cauquenes_sf$wt #pesos x <- temp_cauquenes_sf$t_mean1 # variable de interés global_moran(x, nb, wt)[1] #Extracción I ``` ``` ## $I ## [1] 0.493794 ``` --- ## Estacionariedad y autocorrelación espacial -- - Ahora, para ver el cambio de correlación en el espacio emplearemos la función `hscat()` de gstat para ver si podemos considerar estacionariedad de segundo orden. -- - Esto considerando que al movernos en el espacio la correlación no varÃe significativamente. ```r hscat(formula = t_mean1 ~ 1, #especificamos la variable dependiente data = temp_cauquenes_sf, #tabla donde se resolverá la fórmula breaks = c(1000,2000,5000,7500,10000,15000), #quiebres de distancia ) ``` --- <img src="data:image/png;base64,#DIPGEOPR_02_5_files/figure-html/unnamed-chunk-13-1.png" width="100%" /> --- ## Estacionariedad y autocorrelación espacial -- - Otra forma en la que podemos "evaluar" la estacionariedad de segundo orden es emplear un mapa variogramas. -- - Acá no consideramos la afgrupación de pares de puntos separadamente sino que agrupadamente en una celda determinada ```r vgm_map <- variogram(t_mean1~1, #indicar formula data = temp_cauquenes_sf, #datos donde se resolverá la formula cutoff = 15000, #distancia de corte width = 500, #celda map = T) #indicar salida del mapa #plot del mapa plot(vgm_map, col.regions = viridis::viridis(50), main="Mapa variográfico", xlab="x", ylab="y") ``` --- <img src="data:image/png;base64,#DIPGEOPR_02_5_files/figure-html/unnamed-chunk-15-1.png" width="100%" style="display: block; margin: auto;" /> --- ## BibliografÃa complementaria - Davis, J. (2022). Introduction to Environmental Data Science. Disponible en: https://bookdown.org/igisc/EnvDataSci/ - Fernández, R. & Cotos, T. (2021). EstadÃstica Espacial con R. Disponible en: https://rubenfcasal.github.io/estadistica_espacial/ - Getis, A. (2005). Spatial Pattern Analysis. En: Encyclopedia of Social Measurement, Elsevier, 2005, Pages 627-632. Editor: Kimberly Kempf-Leonard, https://doi.org/10.1016/B0-12-369398-5/00336-4. - Isaaks, E.H., and Srivastava, R.M. (1989). An Introduction to Applied Geostatistics, Oxford University Press, New York. - Olaya, V. (2020). Sistemas de información geográfica. Un libro libre de VÃctor Olaya. Disponible en: https://volaya.github.io/libro-sig/ - Tobler, W.R. (1970). A computer movie simulating urban growth in the Detroit region. Economic Geography, 46, 234-240. DOI: 10.2307/143141. - Webster R., and Oliver M.A. (2001). Geostatistics for Environmental Scientists. John Wiley and Sons, West Sussex, England. --- class: inverse middle 